/**
 * Three circles of equal radius are placed inside a larger circle such that 
 * each pair of circles is tangent to one another and the inner circles do 
 * not overlap. There are four uncovered "gaps" which are to be filled 
 * iteratively with more tangent circles.
 *
 * At each iteration, a maximally sized circle is placed in each gap, which 
 * creates more gaps for the next iteration. After 3 iterations, there are 
 * 108 gaps and the fraction of the area which is not covered by circles is 
 * 0.06790342, rounded to eight decimal places. 
 *
 * What fraction of the area is not covered by circles after 10 iterations?
 * Give your answer rounded to eight decimal places using the format x.xxxxxxxx . 
 *
 * SOLUTION: 
 *
 * This circle packing procedure is known as "Apollonian gasket". The radius
 * of the circles generated in each iteration can be found using Descartes' 
 * theorem.
 *
 * ANSWER: 0.00396087
 */

#include <cmath>
#include <iostream>

static double covered_area(int levels, double b1, double b2, double b3)
{
	// Find a circle tangent to all the three circles. There are in general
	// two solutions. We pick the smaller radius that is positive, i.e. the
	// bigger bend.
	double u = b1 + b2 + b3;
	double v = 2 * sqrt(b1*b2 + b2*b3 + b3*b1);
	double b = u + v;

	// Accumulate the area of this circle and proceed iterations.
	double area = 1.0/(b*b);
	if (levels > 1)
	{
		double a1 = covered_area(levels-1, b, b1, b2);
		double a2 = (b1 == b3)? a1 : covered_area(levels-1, b, b2, b3);
		double a3 = (b2 == b3)? a1 : (b2 == b1)? a2 : covered_area(levels-1, b, b3, b1);
		area += (a1 + a2 + a3);
	}
	return area;
}

void solve_problem_199()
{
	const int levels = 10;

	const double r0 = 1.0;
	const double r1 = (2*sqrt(3.0)-3.0)*r0;

	double area = 3*r1*r1;
	area += covered_area(levels, 1/r1, 1/r1, 1/r1);
	area += 3*covered_area(levels, -1/r0, 1/r1, 1/r1);

	std::cout.precision(8);
	std::cout << std::fixed << (r0*r0 - area) << std::endl;
}
